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The order of the thermal transition in the chiral limit of QCD with two dynamical flavours of 
quarks is a long-standing issue. Still, it is not definitely known whether the transition is of first 
or second order in the continuum limit. Which of the two scenarios is realized has important im¬ 
plications for the QCD phase diagram and the existence of a critical endpoint at finite densities. 
Settling this issue by simulating at successively decreased pion mass was not conclusive yet. Re¬ 
cently, an alternative approach was proposed, extrapolating the first order phase transition found 
at imaginary chemical potential to zero chemical potential with known exponents, which are in¬ 
duced by the Roberge-Weiss symmetry. For staggered fermions on = 4 lattices, this results 
in a first order transition in the chiral limit. Here we report of A^r = 4 simulations with Wilson 
fermions, where the first order region is found to be large. 
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Figure 1: Possible scenarios for the QCD phase diagram at p —0 as function of quark mass. 


1. Introduction 

The clarification of the order of the thermal transition in the chiral limit of Nf = 2 QCD is a 
long-standing issue. Still, it is not definitely known if the transition is of first or second order. This 
is depicted in Figure 1 . Which of the two scenarios is realized has important implications for the 
physical QCD phase diagram, and in particular it is important regarding the possible existence of a 
critical endpoint at finite densities. Settling this issue by simulating at successively decreased pion 
mass was not conclusive yet, primarily because of the increasing demands of the simulations as the 
pion mass is lowered (see e.g. [1] for references and a more detailed introduction to the topic). 

Recently, an alternative approach was proposed [1], which relies on the nontrivial phase struc¬ 
ture of QCD at purely imaginary chemical potential pj. In this region of phase space, one has 
reflection symmetry and extended center (or Roberge-Weiss (RW)) symmetry [4]: 

Z{p) = Zi-p), (1.1) 

Z{p) = Z{p+27iik/N^), k€N . (1.2) 

At imaginary chemical potential p = ipi the sign problem is absent and standard simulation al¬ 
gorithms can be applied. Critical values of Pj = {2k + l)7r/Nc {k G N) mark the boundary be¬ 
tween adjacent center sectors. The transition between these sectors in pi is first order for high and 
crossover for low temperatures. The endpoint of this so-called Roberge-Weiss transition meets with 
the chiraVdeconfinement transition continued from real p. Thus it corresponds to a triple point at 
low and high masses, where first order chiral and deconfinement transitions join it, and to a second 
order endpoint for intermediate masses, where the thermal transition is just a crossover. These 
regions are separated by tricritical points. When mass and Nf are changed at fixed Pj, a phase dia¬ 
gram similar to that shown in Figure 1 for p = 0 emerges. Both are analytically connected when pj 
is varied, which is depicted in Figure 2 (left). More specifically, leaving the critical /i/-values (the 
bottom of the figure), lines of second order transitions depart from the tricritical points, separating 
regions of first order transitions from crossover regions. In the vicinity of the tricritical points, the 
line is governed by tricritical scaling laws, which allows for an extrapolation to the chiral limit. For 
the Nf = 2 backplane, which is of interest in this study, the possible scenarios are shown in Figure 
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Figure 2: Left: Expected QCD phase diagram as function of Right: Possible scenarios in the 

Nf=2 backplane. Both figures follow [1 ]. 


2 (right). If the tricritical point at m = 0 is at negative values of the chiral phase transition is 
second order. On the other hand, if it is at positive values, there exists a first order region at p = 0 
and the transition in the chiral limit must be first order, too. In this way one can clarify the order of 
the chiral limit at zero chemical potential by mapping out the second order line. Using staggered 
fermions on A/j = 4 lattices, it was indeed found that the transition is of first order in the chiral 
limit [1]. These findings have to be contrasted with other fermion discretizations. We report on the 
status of our simulations with Wilson fermions following the same approach. 

2. Simulation details 

We employ the same numerical setup as for the study of the Roberge-Weiss transition de¬ 
scribed in [2], with the standard Wilson gauge action and two flavours of mass-degenerate unim¬ 
proved Wilson fermions. The bare fermion mass m is encapsulated in the hopping parameter 
K = {2(am + 4))^^. Finite temperature on the lattice is given by the inverse temporal extent, 
r = 1/(a'(j3)NT). All simulations were carried out using the OpenCL^ based code CI?QCD [3]^, 
which runs efficiently on Graphic Processing Units (GPUs) on LOEWE-CSC at Goethe-University 
Erankfurt [5] and on E-CSC at GSI Darmstadt [6]. We work at fixed temporal lattice extent = 4, 
leaving the RW-plane Pf = nT investigated in [2]. In the latter study, the lowest mass simulated 
was at fc = 0.165. Here we add tc = 0.170,0.175 and 0.180. In order to locate the critical chemical 
potential for each bare quark mass, simulations at various values of the quark chemical potential 
ap were carried out. Eor each of these parameter sets, temperature scans were carried out with 
Aj3 at least 0.001 around the critical coupling on three spatial volumes, Ncj = 12,16 and 20, corre¬ 
sponding to aspect ratios Ncj/Nt of 3,4 and 5, respectively. After discarding 5k to 10k trajectories 
for thermalization, 40k to 60k trajectories have been simulated on each individual Monte-Carlo 
chain, such that there are at least 100 independent measurements around the critical region. The 

*See www.khronos.org/opencl for more information. 

^Which is now available at github.eom/CL2QCD . 
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Figure 3: Left: Critical coupling as function of p/T for K = 0.170. The line indicates lic{pj). Right: 

Finite size scaling ofB^ and fit for K = 0.165. The vertical lines indicate the fit ranges. 
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Figure 4: Left: Preliminary fit results for K = 0.180. Right: Results for 



autocorrelation on the data was estimated using a python implementation^ of the Wolff method 
[8]. To accumulate statistics faster, we simulated four chains for each parameter set. Observables 
are measured after each trajectory and the acceptance rate in each run was of the order of 75%. 
Additional j3-points have been filled in using Ferrenberg-Swendsen reweighting [9]. 

3. Results 

We define a (pseudo-)critical temperature Tc or coupling jSc by the vanishing of the skewness 

S ={{X-{X))^)/{{X-{X)f)y^ (3.1) 

of a suitable observable X. In this study, we use the chiral condensate (yw) = Nf TrD^^ Figure 3 
(left) shows pc for tc = 0.170. For a fixed K, the PdlL) can be nicely fitted to a quadratic, even func¬ 
tion, as shown in the figure. This allows to interpolate between simulation points and to extrapolate 
towards zero chemical potential. As expected, the results show a decreasing critical temperature as 
the chemical potential approaches zero. The same holds if the (bare) mass is lowered. 

^See github.com/dhesse/py-uwerr . 
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K 

i 8 

a[fm] 

T [MeV] 

mji [MeV] 

0.1800 

4.9519 

0.309(3) 

162(2) 

587(6) 

0.1750 

5.0519 

0.301(3) 

166(2) 

642(7) 

0.1700 

5.1500 

0.288(3) 

174(2) 

699(7) 

0.1650 

5.2420 

0.271(3) 

185(2) 

770(8) 

0.1575 

5.3550 

0.246(3) 

203(2) 

929(10) 


K 

Pc 

an. 

0.165 

5.2421(1) 

0.2408(6) 

0.170 

5.1497(9) 

0.1981(9) 

0.175 

5.0519(3) 

0.1346(17) 

0.180 

4.9520(2) 

0.0090(423) 


Table 1: Results for pc and apf 

from finite size scaling. Table 2: Results of the scale setting. See text for details. 


In the thermodynamic limit V —)■ oo, the Binder cumulant [10] 

B,{X) = {{X - {X)f)/{{X - (X))2)2 (3.2) 

allows to extract the order of the transition. In particular, it takes the values 1 for a first order 
transition and 3 when there is no true phase transition but a crossover. For the case of a second order 
transtion in the 3D Ising universality class, it has a value of around 1.604 (see e.g. [7]). Hence, a 
discontinuity exists when passing from the first order to the crossover region via the second order 
endpoint. To clarify the type of transition on finite lattices one can look at the finite-size scaling of 
the B 4 at pc- In the vicinity of the second order point it scales as [ 10 ] 

B 4 {Pc,N a) = bi+b 2 [{an,)^-{anff]Ny'' + ... . (3.3) 

For the 3D Ising universality class, one has b\ « 1.604 and v « 0.63. The critical coupling n, 
indicates the position of the Z(2) point. 

The values for B 4 {pc,N(y) obtained from the simulations can be fitted to this form. Examples 
are given in Figure 3 (right) and 4 (left). The results for /if are shown in Figure 4 (right). Results 
for pc and /if are summarized in Table 1. As can be seen in the figures, B 4 increases wifh volume 
if fhe transition is a crossover (left of the second order point), whereas it decreases in the first 
order region, ultimately approaching the infinite volume values 3 and 1, respectively. As the (bare) 
quark mass is lowered, /if/T decreases towards zero (see Figure 4 (right)). Note that the results for 
fc = 0.180 are preliminary. With the given data, one has to extrapolate to /if (see Figure 4 (left)). 
Therefore we are cuiTently adding more simulations here in order to get a better estimate of /if. 
Nevertheless, /if is clearly compatible with zero within en'ors given the current data. This indicates 
that also for /i = 0 a first order region is present. 

To relate these findings fo physical scales we performed a series of T = 0 simulafions af or 
close fo fhe respecfive pc for each K. Similar to [2] we generated ^(400) independent configu¬ 
rations on 32^ X 12 lattices. With these at hand, we set the scale by the Wilson flow parameter 
Wo using fhe publicly available code described in [11]. This mefhod is very efficienl and fasf. In 
particular, if is much more precise fhan selling fhe scale by fhe p mass, which was done in [ 2 ]. 
Hence, we revalualed fhe T = 0 simulations from Ibis sludy, carried oul al x: = 0.1575, and in¬ 
clude Ihem here for completeness. In addilion, fhe pion mass was delermined^. The resulls 
for fhe lattice spacing a, fhe critical lemperalure Tc and in physical unils are summarized in 

^We thank F. Depta for carrying out these measurements. 
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Figure 5: Z(2) line in the m,i-{ii/T)^ plane. See text for details. 
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Table 2. The results show that, in terms of pion masses, the first order region is large. Note that 
the lattices coarsen going to lower masses, since j3c decreases. However, all lattices considered are 
very coarse, a > 0.25 fm. Because of this, huge discretization artifacts can be expected. Note that 
mjcL > 5 holds for all our parameter sets, so that finite size effects are negligible. 

4. Discussion & Perspectives 

Our findings are summarized in Figure 5, which shows fhe Z(2) line in fhe plane. 

Leaving fhe RW-plane, the critical line approaches zero jX at clearly non-vanishing pion masses. In 
fact, the lowering of the critical pion masses is steady and relatively mild. This finding corresponds 
to the second scenario sketched in Figure 2 (right). 

For comparison, the results from the previous study using staggered fermions [1] are also 
shown in Figure 5. On coarse lattices, both staggered and Wilson discretizations show similar 
behaviour with a first order region at /x = 0. This region is much larger for Wilson fermions. 
In (^(a)-improved) Wilson studies at N^ = 12 [12] and 16 [13], also shown in the figure, only 
crossover signals are seen, yielding an upper bound for a possible first order region. This suggest 
that the observed wide first order region is to a large extent due to discretization effects. 

To put our N^ = 4 results into perspective, they differ from early results with Wilson fermions 
based on a different simulation strategy [14]. However, they are in accord with modern investiga¬ 
tions. A recent study with ^(a)-improved Wilson-Clover fermions determined a similarly large 
of around 880 MeV for Nf = 3 on N^ = 4 lattices[15]. In another study presented at this conference 
[16], it was shown that the tricritical point in the RW-plane moves to lower masses by only 100 
MeV within our action. Taken together, this may suggest that the ^(a) effects are not yet dominant, 
at least not on N^ = 4 lattices. To reduce cut-off effects, one has to study the NT-dependence of the 
first order region. We will determine this dependence for the jx =0 endpoint of the chiral critical 
line in future studies. 
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